%     VIRTUAL UD - программа для получения эффективных свойств однонаправленного жгута
%     Copyright (C) 2020  Михаил Козлов
%     Сайт проекта: <http://www.fastenerdoctor.ru/>
%     Адрес электронной почты для связи с разработчиком:
%     <doctor.fastener@gmail.com>
% 
%     Эта программа является свободным программным обеспечением: Вы можете
%     распространять ее и/или изменять, соблюдая условия Универсальной общественной
%     лицензии GNU, опубликованной Фондом свободного программного обеспечения,
%     либо редакции 3, либо (на Ваше усмотрение) любой редакции, выпущенной позже.
% 
%     Эта программа распространяется в надежде на то, что она окажется полезной,
%     но БЕЗ КАКИХ-ЛИБО ГАРАНТИЙ, включая даже подразумеваемую гарантию ТОВАРНОЙ ПРИГОДНОСТИ,
%     либо ПРИГОДНОСТИ ДЛЯ ОПРЕДЕЛЕННЫХ ЦЕЛЕЙ. Ознакомьтесь с Универсальной общественной
%     лицензией GNU для получения более подробной информации.
% 
%     Вы должны были получить копию Универсальной общественной лицензии GNU
%     вместе с этой программой. Если Вы ее не получили, то перейдите по адресу:
%     <http://www.gnu.org/licenses/>.
% 
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
% 
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
%     GNU General Public License for more details.
% 
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see <https://www.gnu.org/licenses/>.

function Write_Input_Mesh(Work_Directory, NODES, ELEMENTS, Material_Properties, Coordinate_Systems)

global ErrCode;

node_coord_format = '%.12e';
orientation_format = '%.12e';
mat_prop_format = '%.12e';

fprintf (1, 'Создание .inp файла с сеткой...                         ');

Part_Name = 'UD_RVE';
Job_Name = 'Virtual_UD_RVE_mesh';
Model_Name = 'Virtual_UD_RVE_mesh';

% Проверка нумерации узлов
for i=1:size(NODES,1)
    if (NODES(i,1) ~= i)
        fprintf (1, '\nОшибка функции Write_Input_Mesh: некорректная нумерация узлов!\n');
        ErrCode = -1;
        return;
    end
end

% Проверка нумерации элементов
for i=1:size(ELEMENTS,1)
    if (ELEMENTS(i,1) ~= i)
        fprintf (1, '\nОшибка функции Write_Input_Mesh: некорректная нумерация элементов!\n');
        ErrCode = -1;
        return;
    end
end

File_name = Job_Name;
File_path = [Work_Directory '\' File_name '.inp'];

Fout = fopen(File_path, 'wt');

if (Fout < 0)
    fprintf (1, '\nОшибка функции Write_Input_Mesh: не удалось открыть для записи файл [%s] !\n', File_path);
    ErrCode = -1;
    return;
end

fprintf (Fout, '*Heading\n');

fprintf (Fout, '** Job name: %s Model name: %s\n', Job_Name, Model_Name);
fprintf (Fout, '** Generated by: Virtual_UD\n');

fprintf (Fout, '*Preprint, echo=NO, model=NO, history=NO, contact=NO\n');

fprintf (Fout, '**\n');
fprintf (Fout, '** PARTS\n');
fprintf (Fout, '**\n');

fprintf (Fout, '*Part, name=%s\n', Part_Name);

fprintf (Fout, '*Node\n');

for i=1:size(NODES, 1)
    fprintf (Fout, ['%i, ' node_coord_format ', ' node_coord_format ', ' node_coord_format '\n'], NODES(i,1), NODES(i,2), NODES(i,3), NODES(i,4));
end

C3D8_Set = zeros(size(ELEMENTS,1),1);
C3D20R_Set = zeros(size(ELEMENTS,1),1);

counters = zeros(2,1);

for i=1:size(ELEMENTS,1)
    
    switch ELEMENTS(i,3)

        case 8

            counters(1) = counters(1) + 1;
            C3D8_Set(counters(1)) = i;

        case 20

            counters(2) = counters(2) + 1;
            C3D20R_Set(counters(2)) = i;

    end
    
end

if (counters(1) > 0)
    C3D8_Set = C3D8_Set(1:counters(1));
else
    C3D8_Set = [];
end

if (counters(2) > 0)
    C3D20R_Set = C3D20R_Set(1:counters(2));
else
    C3D20R_Set = [];
end

if (~isempty(C3D8_Set))
    
    if ((Material_Properties.Matrix_nu > 0.49) || (Material_Properties.Fiber_Nu_23 > 0.49))
        fprintf (Fout, '*Element, type=C3D8IH\n');
    else
        fprintf (Fout, '*Element, type=C3D8I\n');
    end
    
    for i=1:size(C3D8_Set,1)
        fprintf (Fout, '%i, %i, %i, %i, %i, %i, %i, %i, %i\n', ELEMENTS(C3D8_Set(i),1), ELEMENTS(C3D8_Set(i),5), ELEMENTS(C3D8_Set(i),6), ELEMENTS(C3D8_Set(i),7), ...
                                                               ELEMENTS(C3D8_Set(i),8), ELEMENTS(C3D8_Set(i),9), ELEMENTS(C3D8_Set(i),10), ELEMENTS(C3D8_Set(i),11), ...
                                                               ELEMENTS(C3D8_Set(i),12));
    end
    
end

if (~isempty(C3D20R_Set))
    
    if ((Material_Properties.Matrix_nu > 0.49) || (Material_Properties.Fiber_Nu_23 > 0.49))
        fprintf (Fout, '*Element, type=C3D20RH\n');
    else
        fprintf (Fout, '*Element, type=C3D20R\n');
    end
    
    for i=1:size(C3D20R_Set,1)
        fprintf (Fout, '%i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i, %i,\n', ELEMENTS(C3D20R_Set(i),1), ELEMENTS(C3D20R_Set(i),5), ELEMENTS(C3D20R_Set(i),6), ELEMENTS(C3D20R_Set(i),7), ...
                                                                                           ELEMENTS(C3D20R_Set(i),8), ELEMENTS(C3D20R_Set(i),9), ELEMENTS(C3D20R_Set(i),10), ELEMENTS(C3D20R_Set(i),11), ...
                                                                                           ELEMENTS(C3D20R_Set(i),12), ELEMENTS(C3D20R_Set(i),13), ELEMENTS(C3D20R_Set(i),14), ELEMENTS(C3D20R_Set(i),15), ...
                                                                                           ELEMENTS(C3D20R_Set(i),16), ELEMENTS(C3D20R_Set(i),17), ELEMENTS(C3D20R_Set(i),18), ELEMENTS(C3D20R_Set(i),19));
        fprintf (Fout, '   %i, %i, %i, %i, %i\n', ELEMENTS(C3D20R_Set(i),20), ELEMENTS(C3D20R_Set(i),21), ELEMENTS(C3D20R_Set(i),22), ELEMENTS(C3D20R_Set(i),23), ELEMENTS(C3D20R_Set(i),24));
    end
    
end

Fiber_Elements = zeros(size(ELEMENTS,1),1);
Matrix_Elements = zeros(size(ELEMENTS,1),1);

counters = zeros(2,1);

for i=1:size(ELEMENTS,1)
    
    switch ELEMENTS(i,4)
        
        case 1
            
            counters(1) = counters(1) + 1;
            Fiber_Elements(counters(1)) = i;
            
        case 2
            
            counters(2) = counters(2) + 1;
            Matrix_Elements(counters(2)) = i;
        
    end
    
end

if (counters(1) > 0)
    Fiber_Elements = Fiber_Elements(1:counters(1));
else
    Fiber_Elements = [];
end

if (counters(2) > 0)
    Matrix_Elements = Matrix_Elements(1:counters(2));
else
    Matrix_Elements = [];
end

fprintf (Fout, '*Nset, nset=ALL_NODES, generate\n');
fprintf (Fout, ' %i, %i, 1\n', 1, NODES(size(NODES,1),1));

fprintf (Fout, '*Elset, elset=ALL_ELEMENTS, generate\n');
fprintf (Fout, ' %i, %i, 1\n', 1, ELEMENTS(size(ELEMENTS,1),1));

if (~isempty(Fiber_Elements))
    fprintf (Fout, '*Elset, elset=Fiber_Elements\n');
    Write_Set(Fout, Fiber_Elements);
end

if (~isempty(Matrix_Elements))
    fprintf (Fout, '*Elset, elset=Matrix_Elements\n');
    Write_Set(Fout, Matrix_Elements);
end

format_line = [orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format ', ' orientation_format '\n'];

fprintf (Fout, '*Orientation, name=Ori_Coord\n');
fprintf (Fout, format_line, Coordinate_Systems{1}.Point_On_Axis_1(1), Coordinate_Systems{1}.Point_On_Axis_1(2), Coordinate_Systems{1}.Point_On_Axis_1(3), ...
                            Coordinate_Systems{1}.Point_On_Plane_12(1), Coordinate_Systems{1}.Point_On_Plane_12(2), Coordinate_Systems{1}.Point_On_Plane_12(3), ...
                            Coordinate_Systems{1}.Origin(1), Coordinate_Systems{1}.Origin(2), Coordinate_Systems{1}.Origin(3));
fprintf (Fout, '3, 0.\n');

fprintf (Fout, '** Section: Fiber_Sectior\n');
fprintf (Fout, '*Solid Section, elset=Fiber_Elements, orientation=Ori_Coord, material=Fiber_Material\n');
fprintf (Fout, ',\n');

fprintf (Fout, '** Section: Matrix_Section\n');
fprintf (Fout, '*Solid Section, elset=Matrix_Elements, material=Matrix_Material\n');
fprintf (Fout, ',\n');

fprintf (Fout, '*End Part\n');

fprintf (Fout, '**\n');
fprintf (Fout, '** ASSEMBLY\n');
fprintf (Fout, '**\n');

fprintf (Fout, '*Assembly, name=Assembly\n');
fprintf (Fout, '**\n');
fprintf (Fout, '*Instance, name=%s-1, part=%s\n', Part_Name, Part_Name);
fprintf (Fout, '*End Instance\n');
fprintf (Fout, '**\n');
fprintf (Fout, '*End Assembly\n');

fprintf (Fout, '**\n');
fprintf (Fout, '** MATERIALS\n');
fprintf (Fout, '**\n');

fprintf (Fout, '*Material, name=Fiber_Material\n');
fprintf (Fout, '*Density\n');
fprintf (Fout, [mat_prop_format ',\n'], Material_Properties.Fiber_rho);
fprintf (Fout, '*Elastic, type=ENGINEERING CONSTANTS\n');
format_line = [mat_prop_format ', ' mat_prop_format ', ' mat_prop_format ', ' mat_prop_format ', ' mat_prop_format ', ' mat_prop_format ', ' mat_prop_format ', ' mat_prop_format '\n'];
fprintf (Fout, format_line, Material_Properties.Fiber_E_1, Material_Properties.Fiber_E_2, Material_Properties.Fiber_E_3, ...
                            Material_Properties.Fiber_Nu_12, Material_Properties.Fiber_Nu_13, Material_Properties.Fiber_Nu_23, ...
                            Material_Properties.Fiber_G_12, Material_Properties.Fiber_G_13);
fprintf (Fout, [mat_prop_format ',\n'], Material_Properties.Fiber_G_23);
fprintf (Fout, '*Expansion, type=ORTHO\n');
fprintf (Fout, [mat_prop_format ', ' mat_prop_format ', ' mat_prop_format '\n'], Material_Properties.Fiber_alpha_11, Material_Properties.Fiber_alpha_22, Material_Properties.Fiber_alpha_33);

fprintf (Fout, '*Material, name=Matrix_Material\n');
fprintf (Fout, '*Density\n');
fprintf (Fout, [mat_prop_format ',\n'], Material_Properties.Matrix_rho);
fprintf (Fout, '*Elastic\n');
fprintf (Fout, [mat_prop_format ', ' mat_prop_format '\n'], Material_Properties.Matrix_E, Material_Properties.Matrix_nu);
fprintf (Fout, '*Expansion\n');
fprintf (Fout, [mat_prop_format ',\n'], Material_Properties.Matrix_alpha);

fclose(Fout);

fprintf (1, 'OK!\n');

% Сообщение об успешном создании файла с результатами
fprintf (1, '\n    Файл [%s] создан успешно!\n', File_path);
